Show the code
source("utils.R")Dependencies
Setup
Show the code
spe <- readRDS("../data/spe.rds")
#define the Z-stacks that you want to compare
zstack_list <- list("-0.09", "0.01", "0.21")
#define the celltype that you want to compare across the stacks - hereby we assume independence across the z-stacks which is an assumption that can be challenged
celltype_ls <- "OD Mature"
selectZstacks <- function(zstack, spe){
sub <- spe[, spe$sample_id == zstack]
pp <- .ppp(sub, marks = "cluster_id")
return(pp)
}
pp_ls <- lapply(zstack_list, selectZstacks, spe)
names(pp_ls) <- zstack_listThe theory of spatial point patterns is discussed in great detail in (Baddeley, Rubak, and Turner 2015). The book has an accompanying package called spatstat which offers great functionality to the theoretical concepts described in the book (Baddeley and Turner 2005). This chapter relies heavily on both publications.
Concepts and Definitions of Point Processes
Point Process
In point pattern analysis we assume that the patterns we observe are a realisation of a stochastic process called a point process. The inferences we make about the point pattern are based on the point process. For example, if a pattern is created by a Poisson point process it will be evenly distributed in the observation window (Baddeley, Rubak, and Turner 2015, 127).
When considering a pattern with \(m\) multiple types, as we do in the (Moffitt et al. 2018) dataset, there are two very closely related concepts. One can view the pattern as a multitype point pattern, where all the points are sampled from the same point process. The other option is to consider the pattern as a multivariate point pattern, where the points come from \(m\) distinct point processes. The difference between these two views is that in the multitype framework we assume the points to stem from the same point process. In the multivariate framework we assume that the types stem from distinct point processes and therefore we can consider dependencies of one type alone. Whether or not the patterns stem from the same point process depends on the biological question. If we analyse two cell types in one slice of a tissue, we should consider them as being sampled from one point process. However, if we consider the distribution of a cell type in two slices of the same tissue we can have grounds to consider them as distinct processes (Baddeley, Rubak, and Turner 2015, 564–65).
Observation Windows
The most common set up in point pattern analysis is what we call window sampling. Instead of observing the entire pattern we observe a subset of this pattern in the so called window. An example could be different small microscopy windows through which a big tissue slice is observed. In this case, it would be wrong to assume the window to be the convex hull around the observed points because they are just a sample of the bigger point pattern (Baddeley, Rubak, and Turner 2015, 143–45).
There is another concept called the small world model. It assumes that points can only be observed in a finite small world and not beyond these boundaries. When thinking of an entire tissue, this is a very common scenario. Cells can only be observed within the tissue and not beyond. In this case, it would be correct to not assume a rectangular observation window but to use methods to estimate an unknown sampling window such as the Ripley-Rasson estimate of a spatial region (Baddeley, Rubak, and Turner 2015, 144–45; Ripley and Rasson 1977).
Show the code
setRiprasWindows <- function(pp){
Window(pp) <- ripras(pp)
return(pp)
}
#the entire point patterns with the ripras windows
pp <- lapply(pp_ls, setRiprasWindows)
separateMarks <- function(pp){
#split the multitype point process into several single type processes
ppls <- split(pp)
return (ppls)
}
#the point patterns separated by their marks
pp_ls <- lapply(pp, separateMarks)Complete Spatial Randomness
Complete spatial randomness (CSR) is often used as the null model for various point patterns. It is the result of a Poisson process. A completely spatial random process is characterised by two properties, homogeneity and independence, as discussed below (Baddeley, Rubak, and Turner 2015, 132).
Homogeneity
“Homogeneity […] means that the expected number of points falling in a region B should be proportional to its area |B|” (Baddeley, Rubak, and Turner 2015, 132) given a proportionality constant \(\lambda\). The constant \(\lambda\) represents the intensity of the process, i.e., the average number of points in a unit area (Baddeley, Rubak, and Turner 2015, 132–33). :
\[ \mathbb{E}[X\cap B] = \lambda |B|. \label{eq:expected_number_points} \]
Independence
Independence implies that in two (non-overlapping) regions \(A\) and \(B\), the number of points \(n(X\cap A)\) and \(n(X\cap B)\) are independent random variables. In other words, the number of points in region \(A\) does not affect the number of points in region \(B\). (Baddeley, Rubak, and Turner 2015, 133).
Inhomogeneous Poisson Process
A Poisson process that is spatially varying in its average density of points is called inhomogeneous. Here, the average density, \(\lambda(u)\), sometimes known as the intensity function (see below), is a function of the spatial location \(u\). In this case, the expected number of points falling into a region \(B\), \(\mu = n(X\cap B)\), is an integration of the intensity function over that region (Baddeley, Rubak, and Turner 2015, 138).
\[ \mu = \int_{B} \lambda(u) du. \label{eq:expected_number_inhomogeneous} \]
Isotropy
A point process is called isotropic, if its statistical properties are invariant to rotations; a CSR process is both stationary and isotropic (Baddeley, Rubak, and Turner 2015, 147).
Stationarity
“A point process is called stationary if, when we view the process through a window W , its statistical properties do not depend on the location of the window in two-dimensional space.” (Baddeley, Rubak, and Turner 2015, 146). This is the case for any homogeneous point process, where the statistical properties of the pattern are unchanged given shifting of the observation window. This means it is stationary in all statistical properties; first-order properties (e.g. intensity) and second-order properties (e.g. correlation) (Baddeley, Rubak, and Turner 2015, 218). Not all metrics assume stationarity in its full sense. Inhomogeneous metrics only assume second-order / correlation stationarity. That means while the intensity function can vary spatially (first-order stationarity is not given), the estimates of correlation functions (e.g. the inhomogeneous K-function) should be the same in parts of the window (Baddeley, Rubak, and Turner 2015, 689 ff.).
Local scaling
If a process is not correlation stationary, so the estimates of the inhomogeneous metric vary between locations, locally-scaled versions of the metric can be applicable. This means in subregions, the process is still stationary and isotropic, but there is a rescaling factor that can vary across the total process (Baddeley, Rubak, and Turner 2015, 246–47).
We can use a permutation test to test the inhomogeneity assumption. In this scenario, we split the patterns into quadrats and compare the estimated functions between the quadrats. It should be noted that this test depends on the arbitrary definition of the quadrats. Given our chosen patterns are not independent but result as marks from an overall point-pattern, the permutation approach is questionable (Baddeley, Rubak, and Turner 2015, 689–93).
Show the code
permutation_test <- function(pp, mark, split, minpoints) {
pp_sel <- subset(pp, marks %in% mark, drop = TRUE)
rho_est <- rhohat(unmark(pp_sel), "x", method="tr")
lambda <- predict(rho_est)
tesselation <- quantess(unmark(pp_sel), "x", 3)
tesselation_split <- nestsplit(pp_sel, tesselation, ny=split)
plot(tesselation_split, main = mark)
tesselation_split$inten <- factor(as.integer(tesselation_split$f1) <= 1, labels=c("Hi","Lo"))
res.scaled <- studpermu.test(tesselation_split, pts ~ inten, summaryfunction=Kscaled,
minpoints = minpoints)
res.inhom <- studpermu.test(tesselation_split, pts ~ inten, summaryfunction=Kinhom,
lambda=lambda, minpoints = minpoints)
#p-value of the local-scaling test
print(paste0(mark,' local scaling test ', res.scaled$p.value))
#p-value of the inhomogeneity test
print(paste0(mark,' inhomogeneity test ', res.inhom$p.value))
}
lapply(c("Microglia", "OD Mature", "Ependymal"), function(x) permutation_test(pp[['0.01']], x, split = 3, minpoints = 10))The p-value of the test for local scaling for microglia cells is \(<0.05\) which indicates that the assumption of local scaling is rejected. Therefore, the distribution of microglia cells is not a scaled version of an overall distribution pattern. The p-value of the test for inhomogeneity for both microglia cells is \(>0.05\) indicating that the assumption of correlation stationarity is not rejected. In this case we can use the inhomogeneous version of the K-function which assumes correlation stationarity.
For ependymal and OD mature cells however, the p-values for both the local scaling test and the inhomogeneity test are \(>0.05\) which means that for this choice of quadrats both the correlation stationarity assumption and the local scaling assumption can’t be rejected.
As the interpretation of the permutation test is highly dependent on the quadrats, the results should be interpreted with care. Both inhomogeneous and locally scaled versions of the summary functions have support and both offer interesting insights into the spatial pattern. Therefore, we will compare all versions and show what the choice of metrics means for their interpretation.
Intensity
Intensity is the expected density of points per unit area. It can be interpreted as the rate of occurrence or the abundance of events. The intensity represents a first order property because it is related to the expected number of points. More formally the average intensity of a point process is defined as:
\[ \bar{\lambda} = \frac{n(x)}{|W|} \label{eq:average_intensity} \]
As this is an average over the entire window, it is a good estimate for a homogeneous point process (Baddeley, Rubak, and Turner 2015, 157–60)
Estimating Intensity
For a homogeneous point process, the intensity can be estimated in a simplistic way: summing the individual intensities of the marks (Baddeley, Rubak, and Turner 2015, 161).
Show the code
intensityPointProcess <- function(pp,mark) if(mark) intensity(pp) else sum(intensity(pp))
intensityPointProcess(pp_ls[['0.01']], mark = FALSE) %>% round(6)[1] 0.001909Otherwise, we can compute the intensity for each mark individually.
Show the code
intensityPointProcess(pp_ls[['0.01']], mark = TRUE) %>% round(8) Ambiguous Astrocyte Endothelial Ependymal Excitatory Inhibitory
0.00024151 0.00020183 0.00014653 0.00008373 0.00036867 0.00061393
Microglia OD Immature OD Mature Pericytes
0.00003031 0.00006249 0.00014278 0.00001750 Kernel Estimation
In kernel estimation, we try to estimate the intensity function \(\lambda(u)\) of the point process. There is a wide variety of kernel estimators (see (Baddeley, Rubak, and Turner 2015, 168)), but a popular choice is the isotropic Gaussian kernel where the standard deviation corresponds to the smoothing bandwidth (Baddeley, Rubak, and Turner 2015, 168).
Show the code
pp_sel <- subset(pp_ls[['0.01']]$`OD Mature`, drop = TRUE)
Dens <- density(pp_sel, sigma = 100)
plot(Dens, main = 'Kernel Density (OD Mature cells)')
Quadrat Counting
In quadrat counting, all points falling into a given quadrat are counted. This gives an overview on the characteristics of the point pattern, such as correlation stationarity (Baddeley, Rubak, and Turner 2015, 163).
Show the code
Q5 <- quadratcount(pp_ls[['0.01']], nx=8, ny=8)
plot(unmark(pp[['0.01']]), main='Unmarked Point Pattern Quadrats')
plot(Q5, col='black', add=TRUE)Under independence assumptions, the quadrat counts can be used for testing homogeneity, i.e., whether the points are distributed evenly across the quadrats (Baddeley, Rubak, and Turner 2015, 164–65).
Show the code
val <- quadrat.test(pp_ls[['0.01']]$`OD Mature`, 5, alternative="regular", method="MonteCarlo")
val
Conditional Monte Carlo test of CSR using quadrat counts
Test statistic: Pearson X2 statistic
data: pp_ls[["0.01"]]$`OD Mature`
X2 = 635.09, p-value = 1
alternative hypothesis: regular
Quadrats: 25 tiles (irregular windows)A p-value of 1 indicates that the null hypothesis of irregularity can not be rejected strongly. Thus, the point pattern of oligodendrocyts is strongly irregular.
Alternatively, we can inspect departures from the hypothesis that points were generated by a Poisson process. We can identify hotspots and coldspots by comparing the standard error of the relrisk function, which computes nonparamatric estimates of the relative risk by kernel smoothing, to the theoretical null distribution of points. The relative risk is the ratio of spatially varying probablilities of different types (Buller 2020).
Show the code
# select marks
selection <- c('OD Mature', 'Ependymal', 'Microglia')
pp_sel <- subset(pp[['0.01']], marks %in% selection, drop = TRUE)
f1 <- pValuesHotspotMarks(pp_sel)
# Plot significant p-values
plot(f1$p, main = "Significant difference\n to Poisson process alpha = 0.05")
Testing for CSR
Whether or not a point process is completely spatially random (CSR) depends on two characteristics: points need to be distributed homogeneously and they have to be independent of each other (see definitions above). There are various ways to test for CSR, here we show the use-case of the clark-evans test (Baddeley, Rubak, and Turner 2015, 165–66).
Show the code
clarkevans.test(pp_ls[['0.01']]$`OD Mature`)
Clark-Evans test
No edge correction
Z-test
data: pp_ls[["0.01"]]$`OD Mature`
R = 0.77286, p-value < 2.2e-16
alternative hypothesis: two-sidedMultitype point process
Show the code
sub <- spe[, spe$sample_id == "0.01"]
(pp <- .ppp(sub, marks = "cluster_id"))Marked planar point pattern: 6111 points
Multitype, with levels =
Ambiguous Astrocyte Endothelial Ependymal Excitatory Inhibitory Microglia OD
Immature OD Mature Pericytes
window: rectangle = [1222.5635, 3012.4248] x [-3993.535, -2202.755] unitsMultitype and Multivariate viewpoint
A pattern with multiple type of points, e.g. cell types, can be seen in different ways. One the one hand, the multitype approach assumes that the points \(x\) were recorded together with with their labels \(m\) and that they were generated at the same time. The locations and labels therefore have a joint distribution \(P(X,M)\). On the other hand one can assume that the pattern with multiple types of points is a combination of several distinct point patterns, one for each type of point. This is the multivariate approach and the different point patterns \(A\) and \(B\) form a joint distribution \(P(A,B)\). To test if the labels depend on the location one can assume the following null hypotheses (Baddeley, Rubak, and Turner 2015, 565–67):
- complete spatial randomness and independence (CSRI): the points are distributed at random; the type of each points is randomly allocated; independence between points of different types; allocation of the types independently of the other points and of its location.
- random labeling: each point is assigned a type at random independently of its location
- independence of components: the points of different types are independent of each other.
Apart from CSRI is is also important for the analysis if we can assume stationarity, i.e. the statistical properties of the point pattern do not change in the window.
For simplicity, we will focus on three cell types of our point pattern: Ependymal, OD Mature and Microglia.
Show the code
marks(pp) <- factor(marks(pp))
selection <- c('OD Mature', 'Ependymal', 'Microglia')
fov_sel <- c('0.01')
pp_sel <- subset(pp, marks %in% selection, drop = TRUE)
spe_sel <- spe[, spe$sample_id == "0.01" & spe$cluster_id %in% selection]We select on fov, which corresponds to one cut in the frontal plane.
Show the code
pp_sel |> as.data.frame() |>
ggplot(aes(x = x, y = y, color = marks)) +
geom_point() +
theme_minimal() +
coord_fixed() +
scale_color_brewer(palette = "Set1")
The summary of pp (point pattern) object returns general properties, plus intensities, combined and per mark type.
Show the code
summary(pp)Marked planar point pattern: 6111 points
Average intensity 0.001906561 points per square unit
Coordinates are given to 4 decimal places
Multitype:
frequency proportion intensity
Ambiguous 773 0.126493200 2.411670e-04
Astrocyte 646 0.105711000 2.015445e-04
Endothelial 469 0.076746850 1.463225e-04
Ependymal 268 0.043855340 8.361288e-05
Excitatory 1180 0.193094400 3.681463e-04
Inhibitory 1965 0.321551300 6.130571e-04
Microglia 97 0.015873020 3.026287e-05
OD Immature 200 0.032727870 6.239767e-05
OD Mature 457 0.074783180 1.425787e-04
Pericytes 56 0.009163803 1.747135e-05
Window: rectangle = [1222.5635, 3012.4248] x [-3993.535, -2202.755] units
(1790 x 1791 units)
Window area = 3205250 square unitsTo get the overall intensity the individual intensities can be summed up. Assuming that the the multitype process is first order stationary (i.e. each sub-process is stationary) the individual intensities sum up to the intensity of the unmarked point process (Baddeley, Rubak, and Turner 2015, 574ff.).
Show the code
sum(intensity(pp)) == intensity(unmark(pp))[1] TRUEThe stationarity assumption is not appropriate in all cases. To assess first-order stationarity visually, we can plot the kernel density estimates per type.
Show the code
ppls <- split(pp_sel) # split by mark
plot(density(ppls))
Ependymal and OD Mature cells are cleary inhomogeneous, while for Microglia cells it is not so clear and we could assume homogeneity, especially because the window is larger than the tissue and there is a tissue border in the bottom middle.
To further inverstiagte the spatial arrangement of the different cell types we can calculate the relative risk, i.e., the probability of observing a given celltype at a given location. It is calculated using the function relrisk. The bandwidth for smoothing is calculated with bw.relrisk and might need to be adjusted (Baddeley, Rubak, and Turner 2015, 577–83).
The relrisk function further gives us the dominant mark for different regions of the tissue of interest. This could be interesting in the annotation of spatial domains. It indicates at each location, which cell type is most likely to occur.
Show the code
rpd <- relrisk(pp_sel, diggle = TRUE)
dom <- im.apply(rpd, which.max)
dom <- eval.im(factor(dom, levels = seq_along(levels(unique(marks(pp_sel)))),
labels = levels(unique(marks(pp_sel)))))
plot(dom,las=2,main="Dominant mark")
Correlation and spacing
Nearest neighbourhood contingency
To further investigate the spatial distribution of the marks we can investigate the nearest neighbourhood of each cell type. One possibility is to work with nearest neighborhood contingency tables developed by (Dixon 2002). The statistical tests are implemented in the R package dixon (de la Cruz 2008).
The measure of segregation \(S\) is defined in (Dixon 2002) as
\[S_{i,j}= \frac{\log[(N_{i,j}/(N_i−N_{i,j})]}{[(N_i−1)/(N−N_i)]}\] where \(N_i\) is the number of individuals \(i\), \(N_{i,j}\) is the number of individuals of type \(i\) with a nearest neighbor of type \(j\), and \(N\) is the total number of individuals.
A value of \(S=0\) is consistent with random labeling. A value larger than 0 indicates that the two types are more segregated than expected by chance, the larger the value the more segregated. Note that segregated means that it is more likely to expect a neigbour of type \(j\) than by chance. In the case that the neigbour is of the same type this is equivalent to “attraction” of the types. On the other hand if \(S<0\) it indicates that type \(j\) is less likely to be a neigbour than by chance. The P-values are calculated using expected numbers of nearest neighbors under the null hypothesis of random labeling using a Monte-Carlo simulation and assumes an asymptotic \(\chi^2\) distribution.
Show the code
out <- dixon(as.data.frame(pp_sel), nsim = 99)1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20,
21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60,
61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80,
81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98,
99.Show the code
out$tablaZ %>%
arrange(desc(abs(`Z `))) %>%
select(-` p-val.Nobs`) From To Obs.Count Exp. Count S Z p-val.Z
1 Ependymal Ependymal 262 87.16 1.96 20.04 0.0000
2 Ependymal OD Mature 3 149.18 -2.04 -16.87 0.0000
3 OD Mature Ependymal 8 149.18 -1.43 -14.26 0.0000
4 OD Mature OD Mature 380 253.83 0.60 11.43 0.0000
5 Ependymal Microglia 3 31.66 -1.07 -5.66 0.0000
6 Microglia Ependymal 9 31.66 -0.68 -4.92 0.0000
7 Microglia OD Mature 67 53.99 0.25 2.60 0.0094
8 Microglia Microglia 21 11.34 0.32 2.50 0.0124
9 OD Mature Microglia 69 53.99 0.12 2.35 0.0190In this table we see that most Ependymal cells are very clustered, while Microglia are more evenly distributed. Further we see that it is less likely to find a Ependymal cells next to a OD mature cells than by chance.
OD Mature cells show this interesting characteristic that they are clustered in some parts of the tissue and more evenly distributed in other parts of the tissue. This characteristic is not visible in the table. The statistic also considers only the nearest neighbour and ignores neighbours that are further away.
Summary functions for pairs of types
Similar to the simple case without marks, it is possible to estimate summary functions. In particular, summary functions between different marks can be calculated. Note that the canonical functions assume that the multi-type process is stationary.
Cross K-function
The cross K-function is a summary function that measures the average number of points of type j within a distance r of a point of type i. The formula is given by:
\[ K(r) = \frac{1}{\lambda_j} \mathbb{E} [t(u,r,X^{j})|u \in X^{i}], \]
where \(X^{i}\) is the point pattern of type \(i\) and \(t(u,r,X^{j})\) is the number of points of type \(j\) in a circle of radius \(r\) around \(u\) (Baddeley, Rubak, and Turner 2015, 594–95). It is important to remember that the homogeneous cross K-function assumes that the multitype process is stationary. If this is not the case, there is a risk in misinterpreting the results. The problem is the confounding between clustering and inhomogeneity, c.f. (Baddeley, Rubak, and Turner 2015, 151–52)
First, we plot an overview over the cross K function for the different types. As we have seen before the assumption of stationarity might be not valid. We will therefore use the inhomogeneous version of the cross K function.
Show the code
resCross <- calcCrossMetricPerFov(
spe,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'Kcross.inhom',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>Show the code
resCross <- subset(resCross, sample_id %in% fov_sel)
plotCrossMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id')[[1]]
The diagonal of the inhomogeneous cross K-function plot shows the K-function for the different marks (indication of Poisson or non-Poisson point processes). Off-diagonal panels give indication of independence of points when the number of points follows the expected K-function but does not imply that the individual marks follow a Poisson process. If the processes of the types are independent, we assume that they are also uncorrelated.
In the example above, assuming that the process is inhomogeneous, the Ependymals cells appear to be regularly spaced, which seems counter intuitive. However, this is the result of the pattern being inhomogeneous with spatially varying intensity. When accounting for this, the pattern is more regular than expected under an inhomogeneous point process. The estimation of the inhomogeneous cross functions is not straightforward and results change based on the estimation of the local intensity and the edge correction, c.f. (Baddeley, Rubak, and Turner 2015, 605).
Let’s focus a bit more on the relationship between Ependymal and the other two cell types. We will also calculate confidence intervals for the different cross K-functions. We have already seen that our dataset most likely does not satisfy the assumption of stationarity. For this reason, we will calculate the inhomogeneous cross K-function. Note that the option to calculate confidence intervals is not yet implemented in spatialFDA.
Remember that the dashed line represents the assumption of a multitype Poisson process. If the line lies above the dotted line, there is indication of attraction while if the line is below the dotted line there is indication of repulsion. In the plot above we can see that there is indication of attraction between Ependymal and OD Mature cells while there is indication of repulsion between Ependymal and Microglia cells.
Cross L-function
Alternatively the L cross function with similar interpretation can be calculated using the Lcross function (Baddeley, Rubak, and Turner 2015, 596ff).
Show the code
resCross <- calcCrossMetricPerFov(
spe,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'Lcross.inhom',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>
<simpleError: Weights in K-function were infinite or NA>Show the code
resCross <- subset(resCross, sample_id %in% fov_sel)
plotCrossMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id')[[1]]
Mark connection function
The mark connection function is the cross pair-correlation function, i.e. the generalization of the pair correlation function to a multitype point processes, divided by the unmarked pair-correlation function. It can be interpreted as the conditional probability that two points a distance \(r\) apart have labels of type \(i\) and of type \(j\), given the presence of those points (Baddeley, Rubak, and Turner 2015, 596–97).
Show the code
plotCrossAll(pp_sel, "markconnect", "iso") +
scale_y_continuous(limits = c(0, 1))Show the code
resCross <- calcCrossMetricPerFov(
spe = spe_sel,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'markconnect',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)
#resCross <- subset(resCross, sample_id %in% fov_sel)
plotCrossMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id')[[1]]
The dashed lines indicate expected values under random labeling. The values measure dependence (or association) between the different labelled points. Positive values indicate that nearby points are more likely to have different types than expected by chance. This positive association between different cell types does not necessarily imply dependence, as it could be influenced by a negative association between cells of the same type, as it it could be the case for the Microglia cells. Furthermore, as the calculation is based on the \(K\) function, the mark connection function assumes homogenity.
Cross F-function (empty space function), cross G-function (Nearest-neighbor function) and cross J-function
The cross F-function is the cumulative distribution function of the distance from a location to the nearest point of the same type. For each type \(i\), it is defined as:
\[F_i(r) = \mathbb{P}\{d(u,X^{i}\leq r\}.\]
The cross G-function is the cumulative distribution function of the distance from a location to the nearest point of another type and is defined as:
\[G_{ij}(r) = \mathbb{P}\{d(x,X^{(j)} \setminus u \leq r \mid X^{(i)} \ \text{has a point at u}).\]
If the points are independent of each other, the G and F function are identical. Both assume that the process is stationary. There are inhomogenous alternatives, in case the intensity is varying. Then we only assume correlation stationarity.
There exists a difference in the interpretation of the theoretical values of the K-cross and the G-cross function. For the K-cross, the theoretical value indicates independence between marks while for the G-cross the theoretical value is consistent with the assumption that the points of type j are Poisson in addition to being independent of the points of type \(i\) (Baddeley, Rubak, and Turner 2015, 597 ff).
The cross J-function is defined as:
\[J_{ij}(r) = \frac{1-G_{ij}(r)}{1-F_{j}(r)}\]
and summarizes the interpoint dependence between type \(i\) and \(j\). Under the hypothesis of independent components, i.e., that the point processes of each type are independent, the G-function is equivalent to the F-function and the J-function is equal to 1 (Baddeley, Rubak, and Turner 2015, 597 ff).
Dot functions
For each K-, G- and J- function, there also exist dot functions, which measure distances from points of one type to points of any type. These functions allow us to measure the dependence of one mark with all other marks at once. For expample, the K-dot function represents the expected number of an other point within distance \(r\) of a typical point of type \(i\) (Baddeley, Rubak, and Turner 2015, 600 ff).
Show the code
plotCrossAll(pp_sel, "Kdot.inhom", "iso")Show the code
resCross <- calcCrossMetricPerFov(
spe,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'Kdot',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)[1] "OD Mature"
<simpleWarning in Kdot(X = structure(list(window = structure(list(type = "rectangle", xrange = c(-3035.012264, -1242.450271), yrange = c(2820.324434, 4609.586219), units = structure(list(singular = "unit", plural = "units", multiplier = 1), class = "unitname")), class = "owin"), n = 330L, x = c(-2948.588648, -2897.57708, -2869.656103, -2891.24141, -2853.006718, -3010.277574, -2971.954101, -2904.677092, -2887.791265, -2851.791437, -2849.074748, -2888.740698, -3001.803085, -2950.421434, -2909.636307, -2891.839974, -2993.495947, -3013.524167, -2979.043681, -2942.557256, -2916.580302, -2942.573152, -2870.153308, -3014.760719, -2965.791655, -2997.982094, -2915.058426, -2971.194414, -2777.493667, -2705.709586, -2700.426346, -2652.577315, -2722.280767, -2696.287446, -2652.228644, -2647.736406, -2696.29761, -2727.428433, -2728.429955, -2686.017141, -2673.563985, -2694.412134, -2668.909929, -2830.438634, -2695.535222, -2745.565611, -2682.19967, -2833.931051, -2762.727439, -2728.982509, -2676.243835, -2774.560203, -2810.140797, -2784.815473, -2707.84147, -2667.451079, -2664.998593, -2750.028718, -2827.729711, -2736.182561, -2696.16138, -2504.201919, -2548.19189, -2624.044011, -2593.025973, -2532.743363, -2488.007598, -2591.692045, -2499.152178, -2496.629962, -2526.047305, -2523.765679, -2461.582417, -2602.590931, -2581.537901, -2516.877432, -2628.696367, -2622.64497, -2622.246798, -2613.625969, -2609.644597, -2608.181092, -2599.903422, -2592.779881, -2583.948804, -2583.539748, -2569.169882, -2521.698539, -2522.350268, -2505.718266, -2503.069553, -2478.307907, -2467.6025, -2466.025129, -2462.221525, -2450.735643, -2484.995191, -2548.471308, -2527.003744, -2624.620707, -2387.26222, -2377.046714, -2408.671885, -2388.248919, -2371.718117, -2400.23381, -2428.068384, -2353.874704, -2381.315461, -2274.867703, -2390.174816, -2226.145481, -2050.582732, -2203.809148, -2133.914786, -2212.342478, -2205.192536, -2205.596189, -2201.984075, -2198.864461, -2188.68064, -2181.834842, -2170.400812, -2169.690486, -2166.254113, -2160.46512, -2154.268849, -2153.373917, -2138.725275, -2133.493169, -2113.480705, -2114.05935, -2100.910119, -2093.056651, -2091.011138, -2081.636, -2080.731905, -2077.984034, -2075.830594, -2072.56278, -2048.586667, -2113.746816, -2111.550558, -2086.316319, -2056.815762, -2053.144007, -2191.755348, -2152.649342, -2144.33852, -2101.243137, -2071.091963, -2193.416918, -1997.291817, -1927.184001, -1894.825623, -1877.58193, -1965.694727, -1887.452165, -1881.639719, -1962.893407, -1865.000144, -1930.817334, -1864.038162, -2015.35057, -1988.84309, -1987.89471, -1978.128122, -1938.633242, -1889.08222, -1870.046028, -2012.379259, -2006.53979, -1964.164552, -1932.265108, -1910.782552, -1905.466485, -1873.323675, -1911.294723, -1724.780357, -1721.007625, -1678.709103, -1678.387408, -1686.465555, -1765.453122, -1708.799512, -1712.155819, -1663.10485, -1674.102591, -1700.476631, -1795.004844, -1697.856433, -1668.806815, -1803.098029, -1792.569437, -1786.286675, -1776.054993, -1755.874682, -1743.991335, -1726.337292, -1696.10723, -1688.327886, -1674.271414, -1649.46646, -1830.313555, -1694.881055, -1685.624107, -1727.650791, -1690.405684, -1812.951339, -1803.781844, -1790.415442, -1784.018361, -1780.136391, -1780.268901, -1775.885117, -1749.007922, -1745.847979, -1743.063749, -1726.660223, -1693.600276, -1689.141444, -1675.739132, -1664.511233, -1659.870135, -1648.112164, -1690.610776, -1793.166854, -1697.788483, -1676.44659, -1650.296953, -1789.7932, -1706.819991, -1662.406809, -1615.650595, -1605.177957, -1602.529167, -1592.058585, -1573.057787, -1569.804767, -1564.214192, -1551.827367, -1511.595991, -1603.712842, -1590.971353, -1546.159478, -1544.921424, -1525.143847, -1487.659978, -1574.386694, -1505.388734, -1615.835468, -1581.025504, -1581.462525, -1576.829273, -1565.313123, -1565.21983, -1611.277761, -1605.427875, -1599.540292, -1579.82756, -1566.308811, -1509.328217, -1581.126602, -1532.077601, -1616.556592, -1529.38126, -1478.284952, -1529.979598, -1517.208185, -1536.40972, -1503.753746, -1632.208867, -1629.323804, -1619.558297, -1556.804559, -1512.797886, -1498.402708, -1470.052083, -1459.861512, -1402.276493, -1396.145301, -1314.961563, -1406.221183, -1272.161465, -1340.028868, -1394.781083, -1357.697515, -1340.224204, -1292.966101, -1291.968755, -1389.667817, -1404.535662, -1340.309407, -1416.89511, -1333.905108, -1316.763197, -1295.85602, -1250.712018, -1292.58811, -1373.066236, -1283.369355, -1271.450575, -1266.981863, -1395.509992, -1359.589533, -1430.799589, -1281.437757, -1344.877437, -1276.351606, -1343.149557, -1419.950009, -2841.209074, -2644.259827, -2643.69724, -2639.769058, -2637.359732, -2635.277076, -2581.88169, -2416.209518, -2214.633261, -2050.116742, -1641.450665, -1834.650755, -1727.594083, -1639.475301, -1640.007241, -1577.979884, -1554.816016, -1439.151853, -1434.588352), y = c(2835.768383, 2875.669847, 2948.229483, 2843.552431, 2855.166397, 3164.341178, 3163.589689, 3206.256229, 3119.59291, 3152.037516, 3079.129589, 3136.041864, 3359.197613, 3275.093089, 3303.010553, 3284.828044, 3226.153905, 3522.723196, 3529.882938, 3600.703083, 3549.700718, 3693.19986, 3676.155387, 4030.151175, 4130.690415, 4313.0572, 4465.901528, 4500.767216, 4440.646041, 4536.866798, 4456.594788, 4545.091597, 4519.838812, 4448.325096, 4603.144283, 4478.351401, 4473.058113, 4597.291869, 4591.912034, 4443.648231, 4606.313836, 4575.335564, 4453.60095, 4078.588959, 4207.565574, 3836.29171, 3679.773345, 3555.017771, 3435.990569, 3498.863633, 3449.436989, 3545.805157, 3244.198876, 3342.16202, 3246.745027, 3402.469467, 3302.340493, 3391.479877, 3203.352345, 3056.727567, 3024.14462, 2891.54081, 2908.026178, 3267.141301, 3278.748286, 3530.591294, 3461.369411, 3726.31803, 3674.495736, 3921.850861, 4373.473106, 4406.869954, 4381.736425, 4235.158056, 4388.333718, 4387.255086, 4420.520191, 4566.340641, 4575.055469, 4469.78641, 4440.323399, 4547.764771, 4508.04966, 4474.9182, 4472.475054, 4597.312847, 4449.857658, 4434.374512, 4504.111755, 4536.325595, 4586.123853, 4464.732421, 4564.902567, 4600.052072, 4512.815431, 4565.600833, 4494.423598, 4441.756616, 4564.44851, 4435.734983, 4509.86677, 4558.178418, 4355.224017, 4357.495516, 4263.549852, 4200.53303, 4105.754506, 4162.514271, 4196.380751, 3826.325709, 3512.282745, 3439.620227, 3662.796772, 4180.678296, 4155.173709, 4339.705654, 4250.440963, 4291.456064, 4276.492436, 4264.710077, 4325.311288, 4322.00383, 4299.989196, 4345.27344, 4334.769063, 4262.122977, 4231.426482, 4333.352374, 4289.159685, 4268.312473, 4292.60158, 4322.26525, 4267.616377, 4264.201252, 4353.534494, 4324.716751, 4284.672823, 4233.788766, 4315.859019, 4291.396787, 4228.866223, 4349.820328, 4235.399399, 4331.398266, 4356.028053, 4266.70718, 4291.559249, 4248.905388, 4232.109364, 4297.501244, 4327.828693, 4339.247729, 4470.757466, 4481.461776, 4536.60647, 4443.20883, 4563.133146, 4511.454778, 4453.436315, 4508.498938, 4588.46565, 4442.970837, 4462.979794, 4302.229293, 4264.182122, 4247.036257, 4278.63612, 4253.038014, 4396.842054, 4291.082779, 4182.015141, 4201.462222, 4205.895258, 4115.78863, 4040.597709, 3769.537019, 2952.121616, 2982.064466, 2968.294008, 2856.142242, 2845.305802, 2956.50451, 2932.266683, 3149.4759, 3055.060419, 3256.318748, 3289.553561, 3281.203782, 3608.616167, 3978.457672, 4089.040147, 4200.543477, 4373.22094, 4340.295201, 4337.945368, 4287.982459, 4320.276062, 4356.900213, 4358.287614, 4352.440736, 4358.53312, 4260.305943, 4372.732153, 4384.356861, 4235.319511, 4329.236487, 4238.609518, 4346.793865, 4470.073025, 4432.745015, 4477.595025, 4457.240233, 4514.699451, 4477.885773, 4453.830857, 4495.541361, 4437.7741, 4582.007973, 4434.900081, 4536.46749, 4518.239437, 4429.185128, 4537.82398, 4485.330097, 4528.331239, 4509.271887, 4514.792449, 4471.367195, 4452.771376, 4545.340296, 4447.432351, 4541.433353, 4450.252368, 4486.47381, 4544.317978, 4573.904384, 4547.850315, 4608.27241, 4537.341946, 4591.899558, 4448.184792, 4555.98108, 4446.71009, 4497.636344, 4505.808874, 4514.294625, 4586.857206, 4597.682458, 4442.396086, 4541.705664, 4583.234706, 4470.846389, 4530.388419, 4456.692326, 4544.021902, 4578.300394, 4314.034808, 4394.325312, 4362.082783, 4341.885911, 4361.914805, 4409.478619, 4193.339064, 4085.75337, 3651.472679, 3445.918756, 3424.976685, 3374.467451, 3399.351568, 3347.002523, 3243.50162, 3189.17197, 3076.263691, 2872.286101, 2894.178291, 2934.592506, 2980.918672, 2905.627809, 3003.508061, 2835.013182, 2854.482881, 2996.671785, 2855.183296, 2967.095426, 2978.217875, 3202.229051, 3155.119982, 3112.528353, 3060.145818, 3020.849459, 3040.192509, 3086.634439, 3123.015574, 3249.282971, 3395.012925, 3351.869929, 3262.627213, 3405.862712, 3237.107675, 3529.319871, 3556.309906, 3583.396021, 3506.982802, 3686.680304, 3923.837035, 3877.275304, 4149.979278, 4353.442328, 4306.131126, 4440.139228, 4565.624204, 3452.244101, 4485.295487, 4559.048569, 4269.735255, 4377.596222, 3404.738948, 4416.37055, 4416.756402, 4208.629329, 4214.600924, 3514.966271, 4248.751093, 4408.22055, 4341.648892, 4568.999682, 4412.646676, 4406.541082, 3248.74563, 4294.518844), markformat = "vector", marks = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L), levels = "OD Mature", class = "factor")), class = "ppp"), r = c(0, 5.05050505050505, 10.1010101010101, 15.1515151515152, 20.2020202020202, 25.2525252525253, 30.3030303030303, 35.3535353535353, 40.4040404040404, 45.4545454545455, 50.5050505050505, 55.5555555555556, 60.6060606060606, 65.6565656565656, 70.7070707070707, 75.7575757575758, 80.8080808080808, 85.8585858585859, 90.9090909090909, 95.959595959596, 101.010101010101, 106.060606060606, 111.111111111111, 116.161616161616, 121.212121212121, 126.262626262626, 131.313131313131, 136.363636363636, 141.414141414141, 146.464646464646, 151.515151515152, 156.565656565657, 161.616161616162, 166.666666666667, 171.717171717172, 176.767676767677, 181.818181818182, 186.868686868687, 191.919191919192, 196.969696969697, 202.020202020202, 207.070707070707, 212.121212121212, 217.171717171717, 222.222222222222, 227.272727272727, 232.323232323232, 237.373737373737, 242.424242424242, 247.474747474747, 252.525252525253, 257.575757575758, 262.626262626263, 267.676767676768, 272.727272727273, 277.777777777778, 282.828282828283, 287.878787878788, 292.929292929293, 297.979797979798, 303.030303030303, 308.080808080808, 313.131313131313, 318.181818181818, 323.232323232323, 328.282828282828, 333.333333333333, 338.383838383838, 343.434343434343, 348.484848484848, 353.535353535354, 358.585858585859, 363.636363636364, 368.686868686869, 373.737373737374, 378.787878787879, 383.838383838384, 388.888888888889, 393.939393939394, 398.989898989899, 404.040404040404, 409.090909090909, 414.141414141414, 419.191919191919, 424.242424242424, 429.292929292929, 434.343434343434, 439.393939393939, 444.444444444444, 449.49494949495, 454.545454545455, 459.59595959596, 464.646464646465, 469.69696969697, 474.747474747475, 479.79797979798, 484.848484848485, 489.89898989899, 494.949494949495, 500)): input string '<b7>' cannot be translated to UTF-8, is it valid in 'ANSI_X3.4-1968'?>
[1] "Ependymal"
[1] "Microglia"Show the code
resCross <- subset(resCross, sample_id %in% fov_sel)
plotMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id', ID = "ID")
The dot functions are useful summary statistics to analyse the dependence of one mark with all other marks.
Summary function within and between types
In our original dataset, we have a large number of different marks. We picked three: OD mature, Ependymal and Microglia for illustrative purposes. An alternative to looking at all cross summary function combinations, it is possible to compare between and within types (Baddeley, Rubak, and Turner 2015).
Mark equality function
The Mark or Type Equality function for a stationary multitype point process measures the correlation between types of two points separated by distance r. It is the sum of the mark connection function of all pairs of points of the same type.
If k < 1, points at distance r are less likely than expected to be of the same type. If > 1, they are more likely to be of the same type. The value 1 indicates a lack of correlation (Baddeley, Rubak, and Turner 2015, 603 ff).
Show the code
resCross <- calcMetricPerFov(
spe,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'markcorr',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)
resCross <- subset(resCross, sample_id %in% fov_sel)
plotMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id')
We can see that in our dataset that it the more likely it is to find points of the same type at shorter distances. The curve never crosses the dashed line at 1, which means that it is generally more likely to find points of the same type at any distance than expected by chance.
Tests of randomness and independence
In a multitype point process, there are usually two interesting hypotheses:
- random-labeling hypothesis: the allocation of labels to the points is random
- independent component hypothesis: there is independence between different type of points
If both statments are correct, the point pattern is considered to be complete spatially random and independent (CSRI), the marked analog to complete spatial randomness (CSR) (Baddeley, Rubak, and Turner 2015, 605 ff).
Testing random labelling
The random labeling test is most logical when the marks represents its status, which is not most appropriate assumption when considering cell types. Testing for random labeling can be done using permutation test, in which the labels are randomly permuted. Random labeling can be assumed if the permuted datasets are statistically equivalent to the original dataset (Baddeley, Rubak, and Turner 2015, 609 ff).
Testing the indepenence of components assumption
The i-to-j functions are useful to test the independence of different subprocesses. If the processes of type i and j are independent, then \(K_{ij} = \pi r^2, G_{ij}(r) = F_{j}(r), J_{ij}(r) \equiv 1\). Alternatively, randomization tests can be used in which simulated patterns from the dataset are generated and randomly split into subpatterns. These are then compared to the null hypothesis in which all subpatterns should be statistically equivalent to the original. However, this approach assumes stationarity and there is a need to handle edge effects (Baddeley, Rubak, and Turner 2015, 606 ff).
Assuming stationarity of the total pattern
As outlined above, the homogeneous cross correlation and spacing functions assume stationarity, whereas the inhomogeneous functions only assume correlation stationarity. First-order stationarity is not given in our dataset, when we look at the different patterns individually. However, using the total (unmarked) pattern, we could assume first-order stationarity, since the intensity is the same across the pattern.
Show the code
plot(density(unmark(pp)))
Let’s look at the homogeneous cross K-function.
Show the code
resCross <- calcCrossMetricPerFov(
spe,
selection = c("OD Mature", "Ependymal", "Microglia"),
subsetby = 'sample_id',
fun = 'Kcross',
marks = 'cluster_id',
r_seq = seq(0, 500, length.out = 100),
by = c('Animal_ID', 'sample_id')
)
resCross <- subset(resCross, sample_id %in% fov_sel)
plotCrossMetricPerFov(resCross, theo = TRUE, correction = "iso", x = "r", image_id = 'sample_id')[[1]]
The result is different from the previous analysis. The Ependymal cells now appear to be clustered. This is because stationarity assumes that Ependymal cells could theoretically be present in the total observation window. If this assumption is justified, depends on the context and research question. The interpretation of the results should always be done with the assumption of stationarity or inhomogeneity in mind and should be reported in an analysis.
Session info
Show the code
sessionInfo()R version 4.3.1 (2023-06-16)
Platform: aarch64-apple-darwin20 (64-bit)
Running under: macOS Sonoma 14.5
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.11.0
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
time zone: Europe/Zurich
tzcode source: internal
attached base packages:
[1] stats4 stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] magrittr_2.0.3 stringr_1.5.0
[3] dixon_0.0-8 splancs_2.01-44
[5] spdep_1.2-8 spData_2.3.0
[7] tmap_3.3-4 scater_1.28.0
[9] scran_1.28.2 scuttle_1.10.3
[11] SFEData_1.2.0 SpatialFeatureExperiment_1.2.3
[13] Voyager_1.2.7 rgeoda_0.0.10-4
[15] digest_0.6.33 ncf_1.3-2
[17] sf_1.0-16 reshape2_1.4.4
[19] patchwork_1.2.0 STexampleData_1.8.0
[21] ExperimentHub_2.8.1 AnnotationHub_3.8.0
[23] BiocFileCache_2.8.0 dbplyr_2.3.4
[25] RANN_2.6.1 seg_0.5-7
[27] sp_2.1-1 rlang_1.1.1
[29] ggplot2_3.5.1 dplyr_1.1.3
[31] mixR_0.2.0 spatstat_3.0-6
[33] spatstat.linnet_3.1-1 spatstat.model_3.2-6
[35] rpart_4.1.19 spatstat.explore_3.2-3
[37] nlme_3.1-162 spatstat.random_3.1-6
[39] spatstat.geom_3.2-5 spatstat.data_3.0-1
[41] SpatialExperiment_1.10.0 SingleCellExperiment_1.22.0
[43] SummarizedExperiment_1.30.2 Biobase_2.60.0
[45] GenomicRanges_1.52.1 GenomeInfoDb_1.36.4
[47] IRanges_2.34.1 S4Vectors_0.38.2
[49] BiocGenerics_0.46.0 MatrixGenerics_1.12.3
[51] matrixStats_1.0.0 spatialFDA_0.99.0
loaded via a namespace (and not attached):
[1] spatstat.sparse_3.0-2 bitops_1.0-7
[3] httr_1.4.7 RColorBrewer_1.1-3
[5] tools_4.3.1 utf8_1.2.3
[7] R6_2.5.1 HDF5Array_1.28.1
[9] mgcv_1.9-1 rhdf5filters_1.12.1
[11] withr_2.5.1 gridExtra_2.3
[13] leaflet_2.2.0 leafem_0.2.3
[15] cli_3.6.1 labeling_0.4.3
[17] proxy_0.4-27 R.utils_2.12.2
[19] dichromat_2.0-0.1 scico_1.5.0
[21] limma_3.56.2 rstudioapi_0.15.0
[23] RSQLite_2.3.1 generics_0.1.3
[25] crosstalk_1.2.0 Matrix_1.5-4.1
[27] ggbeeswarm_0.7.2 fansi_1.0.5
[29] abind_1.4-5 R.methodsS3_1.8.2
[31] terra_1.7-55 lifecycle_1.0.3
[33] yaml_2.3.7 edgeR_3.42.4
[35] rhdf5_2.44.0 tmaptools_3.1-1
[37] grid_4.3.1 blob_1.2.4
[39] promises_1.2.1 dqrng_0.3.1
[41] crayon_1.5.2 lattice_0.21-8
[43] beachmat_2.16.0 KEGGREST_1.40.1
[45] magick_2.8.0 pillar_1.9.0
[47] knitr_1.44 metapod_1.7.0
[49] rjson_0.2.21 boot_1.3-28.1
[51] codetools_0.2-19 wk_0.8.0
[53] glue_1.6.2 vctrs_0.6.4
[55] png_0.1-8 gtable_0.3.4
[57] cachem_1.0.8 xfun_0.40
[59] S4Arrays_1.0.6 mime_0.12
[61] DropletUtils_1.20.0 units_0.8-4
[63] statmod_1.5.0 bluster_1.10.0
[65] interactiveDisplayBase_1.38.0 ellipsis_0.3.2
[67] bit64_4.0.5 filelock_1.0.2
[69] irlba_2.3.5.1 vipor_0.4.5
[71] KernSmooth_2.23-21 colorspace_2.1-0
[73] DBI_1.1.3 raster_3.6-26
[75] tidyselect_1.2.0 bit_4.0.5
[77] compiler_4.3.1 curl_5.1.0
[79] BiocNeighbors_1.18.0 DelayedArray_0.26.7
[81] scales_1.3.0 classInt_0.4-10
[83] rappdirs_0.3.3 goftest_1.2-3
[85] fftwtools_0.9-11 spatstat.utils_3.0-5
[87] rmarkdown_2.25 XVector_0.40.0
[89] base64enc_0.1-3 htmltools_0.5.6.1
[91] pkgconfig_2.0.3 sparseMatrixStats_1.12.2
[93] fastmap_1.1.1 htmlwidgets_1.6.2
[95] shiny_1.7.5.1 DelayedMatrixStats_1.22.6
[97] farver_2.1.1 jsonlite_1.8.7
[99] BiocParallel_1.34.2 R.oo_1.25.0
[101] BiocSingular_1.16.0 RCurl_1.98-1.12
[103] GenomeInfoDbData_1.2.10 s2_1.1.4
[105] Rhdf5lib_1.22.1 munsell_0.5.0
[107] Rcpp_1.0.11 ggnewscale_0.4.9
[109] viridis_0.6.4 stringi_1.7.12
[111] leafsync_0.1.0 zlibbioc_1.46.0
[113] plyr_1.8.9 parallel_4.3.1
[115] ggrepel_0.9.4 deldir_1.0-9
[117] Biostrings_2.68.1 stars_0.6-4
[119] splines_4.3.1 tensor_1.5
[121] locfit_1.5-9.8 igraph_1.5.1
[123] ScaledMatrix_1.8.1 BiocVersion_3.17.1
[125] XML_3.99-0.14 evaluate_0.22
[127] BiocManager_1.30.22 httpuv_1.6.11
[129] tidyr_1.3.0 purrr_1.0.2
[131] polyclip_1.10-6 rsvd_1.0.5
[133] lwgeom_0.2-13 xtable_1.8-4
[135] e1071_1.7-13 RSpectra_0.16-1
[137] later_1.3.1 viridisLite_0.4.2
[139] class_7.3-22 tibble_3.2.1
[141] memoise_2.0.1 beeswarm_0.4.0
[143] AnnotationDbi_1.62.2 cluster_2.1.4 ©2024 The pasta authors. Content is published under Creative Commons CC-BY-4.0 License for the text and GPL-3 License for any code.
References
Baddeley, Adrian, Ege Rubak, and Rolf Turner. 2015. Spatial Point Patterns. 1st ed. CRC Interdisciplinary Statistics Series. CRC Press, Taylor & Francis Group.
Baddeley, Adrian, and Rolf Turner. 2005. “Spatstat: An R Package for Analyzing Spatial Point Patterns.” Journal of Statistical Software 12 (January): 1–42. https://doi.org/10.18637/jss.v012.i06.
Buller. 2020. “Areas of a Spatial Segregation Model Significantly Different from Null Expectations.” Ian Buller, Ph.D., M.A. https://idblr.rbind.io/post/pvalues-spatial-segregation/.
de la Cruz, Marcelino. 2008. “Métodos Para Analizar Datos Puntuales.” In Introduccion Al Analisis Espacial de Datos En Ecologia y Ciencias Ambientales: Metodos y Aplicaciones, 76–127.
Dixon, Philip M. 2002. “Nearest-Neighbor Contingency Table Analysis of Spatial Segregation for Several Species.” Écoscience 9 (2): 142–51. https://doi.org/10.1080/11956860.2002.11682700.
Moffitt, Jeffrey R., Dhananjay Bambah-Mukku, Stephen W. Eichhorn, Eric Vaughn, Karthik Shekhar, Julio D. Perez, Nimrod D. Rubinstein, et al. 2018. “Molecular, Spatial, and Functional Single-Cell Profiling of the Hypothalamic Preoptic Region.” Science 362 (6416): eaau5324. https://doi.org/10.1126/science.aau5324.
Ripley, B. D., and J.-P. Rasson. 1977. “Finding the Edge of a Poisson Forest.” Journal of Applied Probability 14 (3): 483–91. https://doi.org/10.2307/3213451.